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Continuous-Time Signals 


Signals occur in a wide range of physical phenomenon. They might be human speech, 
blood pressure variations with time, seismic waves, radar and sonar signals, pictures or 
images, stress and strain signals in a building structure, stock market prices, a city's 
population, or temperature across a plate. These signals are often modeled or represented 
by areal or complex valued mathematical function of one or more variables. For example, 
speech is modeled by a function representing air pressure varying with time. The function 
is acting as a mathematical analogy to the speech signal and, therefore, is called an analog 
signal. For these signals, the independent variable is time and it changes continuously so 
that the term continuous-time signal is also used. In our discussion, we talk of the 
mathematical function as the signal even though it is really a model or representation of 
the physical signal. 


The description of signals in terms of their sinusoidal frequency content has proven to be 
one of the most powerful tools of continuous and discrete-time signal description, 
analysis, and processing. For that reason, we will start the discussion of signals with a 
development of Fourier transform methods. We will first review the continuous-time 
methods of the Fourier series (FS), the Fourier transform or integral (FT), and the Laplace 
transform (LT). Next the discrete-time methods will be developed in more detail with the 
discrete Fourier transform (DFT) applied to finite length signals followed by the discrete- 
time Fourier transform (DTFT) for infinitely long signals and ending with the Z-transform 
which allows the powerful tools of complex variable theory to be applied. 


More recently, a new tool has been developed for the analysis of signals. Wavelets and 
wavelet transforms [link], [link], [link], [link], [link] are another more flexible expansion 
system that also can describe continuous and discrete-time, finite or infinite duration 
signals. We will very briefly introduce the ideas behind wavelet-based signal analysis. 


The Fourier Series 


The problem of expanding a finite length signal in a trigonometric series was posed and 
studied in the late 1700's by renowned mathematicians such as Bernoulli, d'Alembert, 
Euler, Lagrange, and Gauss. Indeed, what we now call the Fourier series and the formulas 
for the coefficients were used by Euler in 1780. However, it was the presentation in 1807 
and the paper in 1822 by Fourier stating that an arbitrary function could be represented by 
a series of sines and cosines that brought the problem to everyone's attention and started 
serious theoretical investigations and practical applications that continue to this day [link], 
[link], [link], [link], [link], [link]. The theoretical work has been at the center of analysis 
and the practical applications have been of major significance in virtually every field of 
quantitative science and technology. For these reasons and others, the Fourier series is 
worth our serious attention in a study of signal processing. 


Definition of the Fourier Series 


We assume that the signal x(t) to be analyzed is well described by a real or complex 
valued function of a real variable t defined over a finite interval {0 < t < T}. The 
trigonometric series expansion of x(t) is given by 

Equation: 
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where xx (t) = cos (2akt/T) and y;, (t) = sin (2akt/T) are the basis functions for the 
expansion. The energy or power in an electrical, mechanical, etc. system is a function of 
the square of voltage, current, velocity, pressure, etc. For this reason, the natural setting for 
a representation of signals is the Hilbert space of L? [0, T']. This modern formulation of 
the problem is developed in [link], [link]. The sinusoidal basis functions in the 
trigonometric expansion form a complete orthogonal set in L? (0, T|. The orthogonality is 
easily seen from inner products 

Equation: 
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and 
Equation: 
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where 6(t) is the Kronecker delta function with 6(0) = 1 and 6(k 4 0) = O. Because of 
this, the kth coefficients in the series can be found by taking the inner product of x(t) with 
the kth basis functions. This gives for the coefficients 

Equation: 


and 
Equation: 
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where 7’ is the time interval of interest or the period of a periodic signal. Because of the 
orthogonality of the basis functions, a finite Fourier series formed by truncating the 
infinite series is an optimal least squared error approximation to x(t). If the finite series is 
defined by 

Equation: 
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the squared error is 
Equation: 
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which is minimized over all a(k) and b(k) by [link] and [link]. This is an extraordinarily 
important property. 


It follows that if x (t) € L? [0,7], then the series converges to x(t) in the sense that 

€ — Oas N — oo[link], [link]. The question of point-wise convergence is more difficult. 
A sufficient condition that is adequate for most application states: If f(a) is bounded, is 
piece-wise continuous, and has no more than a finite number of maxima over an interval, 
the Fourier series converges point-wise to f(x) at all points of continuity and to the 
arithmetic mean at points of discontinuities. If f(x) is continuous, the series converges 
uniformly at all points [link], [link], [link]. 


A useful condition [link], [link] states that if x(t) and its derivatives through the gth 
derivative are defined and have bounded variation, the Fourier coefficients a(k) and b(k) 
asymptotically drop off at least as fast as car as k — oo. This ties global rates of 
convergence of the coefficients to local smoothness conditions of the function. 


The form of the Fourier series using both sines and cosines makes determination of the 
peak value or of the location of a particular frequency term difficult. A different form that 
explicitly gives the peak value of the sinusoid of that frequency and the location or phase 
shift of that sinusoid is given by 

Equation: 
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and, using Euler's relation and the usual electrical engineering notation of 7 = ./—1, 
Equation: 
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the complex exponential form is obtained as 


Equation: 
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where 
Equation: 


The coefficient equation is 
Equation: 


The coefficients in these three forms are related by 
Equation: 
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and 
Equation: 
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It is easier to evaluate a signal in terms of c(k) or d(k) and 0(k) than in terms of a(k) and 
b(k). The first two are polar representation of a complex value and the last is rectangular. 


The exponential form is easier to work with mathematically. 


Although the function to be expanded is defined only over a specific finite region, the 
series converges to a function that is defined over the real line and is periodic. It is equal 
to the original function over the region of definition and is a periodic extension outside of 
the region. Indeed, one could artificially extend the given function at the outset and then 
the expansion would converge everywhere. 


A Geometric View 


It can be very helpful to develop a geometric view of the Fourier series where x(t) is 
considered to be a vector and the basis functions are the coordinate or basis vectors. The 
coefficients become the projections of x(t) on the coordinates. The ideas of a measure of 
distance, size, and orthogonality are important and the definition of error is easy to picture. 
This is done in [link], [link], [link] using Hilbert space methods. 


Properties of the Fourier Series 


The properties of the Fourier series are important in applying it to signal analysis and to 
interpreting it. The main properties are given here using the notation that the Fourier series 
of a real valued function x(t) over {0 < t < T} is given by F{x(t)} = c(k) and Z (t) 
denotes the periodic extensions of x(t). 


1. Linear: F{x + y} = F{x} + F{y} 

Idea of superposition. Also scalability: F{ar} = aF {x} 
2. Extensions of x(t): ¥(t) = %(t + T) 

£ (t) is periodic. 
3. Even and Odd Parts: x(t) = u(t) + ju(t) and 

C (k) = A(k) + jB(k) = |C(k)| e™) 


U v A B |C| 6 
even 0 even 0 even 0 


odd 0 0 odd even 0 


0 even 0 even even 1/2 


0 odd odd 0 even m/2 


4. Convolution: If continuous cyclic convolution is defined by 
Equation: 
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then F{h(t) oa(t)} = F{h(t)} F{x(t)} 
5. Multiplication: If discrete convolution is defined by 
Equation: 


then A{h(t) x(t)} = F{h(t)}*F{ax(t)} 
This property is the inverse of property 4 and vice versa. 
.Parseval: + fy |a (é)|?dt = Te_. |C (A) 
This property says the energy calculated in the time domain is the same as that 
calculated in the frequency (or Fourier) domain. 
. Shift: F {% (t — to)} = C(k) e— Patok/T 
A shift in the time domain results in a linear phase shift in the frequency domain. 
8. Modulate: F {x (t)e?"™*/T} = C(k—K) 
Modulation in the time domain results in a shift in the frequency domain. This 
property is the inverse of property 7. 
9. Orthogonality of basis functions: 
Equation: 
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Orthogonality allows the calculation of coefficients using inner products in [link] and 


[link]. It also allows Parseval's Theorem in property 6. A relaxed version of 
orthogonality is called “tight frames" and is important in over-specified systems, 
especially in wavelets. 


Examples 


e An example of the Fourier series is the expansion of a square wave signal with period 
27. The expansion is 
Equation: 
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Because x(t) is odd, there are no cosine terms (all a(k) = 0) and, because of its 
symmetries, there are no even harmonics (even k terms are zero). The function is 


well defined and bounded; its derivative is not, therefore, the coefficients drop off as 
1 


k . 
e A second example is a triangle wave of period 27. This is a continuous function 
where the square wave was not. The expansion of the triangle wave is 


Equation: 
; I 3 1. 
z(t) = —|sin (t) — az ain (3t) + m2 sin (5t) +--+]. 


Here the coefficients drop off as er since the function and its first derivative exist 
and are bounded. 


Note the derivative of a triangle wave is a square wave. Examine the series coefficients to 
see this. There are many books and web sites on the Fourier series that give insight 
through examples and demos. 


Theorems on the Fourier Series 


Four of the most important theorems in the theory of Fourier analysis are the inversion 
theorem, the convolution theorem, the differentiation theorem, and Parseval's theorem 
[link]. 


e The inversion theorem is the truth of the transform pair given in [link], [link], and 
[link].. 

e The convolution theorem is property 4. 

e The differentiation theorem says that the transform of the derivative of a function is 
jw times the transform of the function. 

e Parseval's theorem is given in property 6. 


All of these are based on the orthogonality of the basis function of the Fourier series and 
integral and all require knowledge of the convergence of the sums and integrals. The 
practical and theoretical use of Fourier analysis is greatly expanded if use is made of 
distributions or generalized functions (e.g. Dirac delta functions, 5(¢)) [link], [link]. 
Because energy is an important measure of a function in signal processing applications, 


the Hilbert space of L? functions is a proper setting for the basic theory and a geometric 
view can be especially useful [link], [link]. 


The following theorems and results concern the existence and convergence of the Fourier 
series and the discrete-time Fourier transform [link]. Details, discussions and proofs can 
be found in the cited references. 


If f(z) has bounded variation in the interval (—7, 7), the Fourier series 
corresponding to f(x) converges to the value f(x) at any point within the interval, at 
which the function is continuous; it converges to the value [f(z +0) + f (x — 0)] 
at any such point at which the function is discontinuous. At the points 7, —7 it 
converges to the value >[f (—m + 0) + f (m — 0)]. [link] 

If f(x) is of bounded variation in (—7, 7), the Fourier series converges to f(z), 
uniformly in any interval (a, b) in which f(z) is continuous, the continuity at a and b 
being on both sides. [link] 

If f(z) is of bounded variation in (—7, 7), the Fourier series converges to 

+(f (« + 0) + f (a — 0)], bounded throughout the interval (—7, 7). [link] 

If f(z) is bounded and if it is continuous in its domain at every point, with the 
exception of a finite number of points at which it may have ordinary discontinuities, 
and if the domain may be divided into a finite number of parts, such that in any one 
of them the function is monotone; or, in other words, the function has only a finite 
number of maxima and minima in its domain, the Fourier series of f(x) converges to 
f(a) at points of continuity and to +[f (a + 0) + f (x — 0)} at points of 
discontinuity. [link], [link] 

If f(a) is such that, when the arbitrarily small neighborhoods of a finite number of 
points in whose neighborhood | f(x)| has no upper bound have been excluded, f(x) 
becomes a function with bounded variation, then the Fourier series converges to the 
value +[f (x + 0) + f (x — 0)], at every point in (—7, 7), except the points of 
infinite discontinuity of the function, provided the improper integral f"_ f (x)dzx 
exist, and is absolutely convergent. [link] 

If f is of bounded variation, the Fourier series of f converges at every point x to the 
value [f(a + 0) + f(x — 0)]/2. If f is, in addition, continuous at every point of an 
interval I = (a, b), its Fourier series is uniformly convergent in I. [link] 

If a(k) and b(k) are absolutely summable, the Fourier series converges uniformly to 
f(a) which is continuous. [link] 

If a(k) and b(k) are square summable, the Fourier series converges to f(x) where it 
is continuous, but not necessarily uniformly. [link] 

Suppose that f(x) is periodic, of period X, is defined and bounded on [0, X] and that 
at least one of the following four conditions is satisfied: (i) f is piecewise monotonic 
on [0, X], (ii) f has a finite number of maxima and minima on [0, X] and a finite 
number of discontinuities on [0, X], (iii) f is of bounded variation on |0, X], (iv) f is 
piecewise smooth on [0, X]: then it will follow that the Fourier series coefficients 


may be defined through the defining integral, using proper Riemann integrals, and 
that the Fourier series converges to f(z) at a.a.x, to f(x) at each point of continuity 
of f, and to the value +[f (x~) + f (x*)] at all «. [Link] 

¢ For any 1 < p< co andany f € C” (8"), the partial sums 
Equation: 
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converge to f, uniformly as n — oo; in fact, ||; — f||,,, is bounded by a constant 
multiple of n-?+1/?. [link] 


The Fourier series expansion results in transforming a periodic, continuous time function, 
£ (t), to two discrete indexed frequency functions, a(k) and b(k) that are not periodic. 


The Fourier Transform 


Many practical problems in signal analysis involve either infinitely long or very long 
signals where the Fourier series is not appropriate. For these cases, the Fourier transform 
(FT) and its inverse (IFT) have been developed. This transform has been used with great 
success in virtually all quantitative areas of science and technology where the concept of 
frequency is important. While the Fourier series was used before Fourier worked on it, the 
Fourier transform seems to be his original idea. It can be derived as an extension of the 
Fourier series by letting the length or period 7’ increase to infinity or the Fourier transform 
can be independently defined and then the Fourier series shown to be a special case of it. 
The latter approach is the more general of the two, but the former is more intuitive [Link], 
[link]. 


Definition of the Fourier Transform 


The Fourier transform (FT) of a real-valued (or complex) function of the real-variable ¢ is 
defined by 
Equation: 


X (w) = [- a(t)e % dt 


co 


giving a complex valued function of the real variable w representing frequency. The 
inverse Fourier transform (IFT) is given by 
Equation: 


1 - 
at) = | X (w) e%* du. 


Because of the infinite limits on both integrals, the question of convergence is important. 
There are useful practical signals that do not have Fourier transforms if only classical 
functions are allowed because of problems with convergence. The use of delta functions 
(distributions) in both the time and frequency domains allows a much larger class of 
signals to be represented [link]. 


Properties of the Fourier Transform 


The properties of the Fourier transform are somewhat parallel to those of the Fourier series 
and are important in applying it to signal analysis and interpreting it. The main properties 
are given here using the notation that the FT of a real valued function x(t) over all time ¢ 
is given by F{x} = X(w). 


1. Linear: F{x + y} = F{x} + F{y} 
2. Even and Oddness: if 2(t) = u(t) + ju(t) and X(w) = A(w) + jB(w) then 


U v A B |X| 6 
even 0 even 0 even 0 
odd 0 0 odd even 0 
0 even 0 even even m/2 
0 odd odd 0 even 1/2 


3. Convolution: If continuous convolution is defined by: 
Equation: 
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then F{h(t)*a(t)} = F{h(t)} F{2(t)} 


4. Multiplication: es att = = a )}*F {ax (t)} 
5. Parseval: [°° |x (t)|?dt = bt f° |X (w)|?dw 
6. Shift: F {x na —T)}= X(w He oF 
7. Modulate: F {x (t)e?"*') = X (w — 2K) 
8. Derivative: F {2} = jwX (w) 
9. Stretch: F {x (at)} = Taj X (w/a) 
10. Orthogonality: fees e Iitejwat — Id (wy — wa) 


Examples of the Fourier Transform 


Deriving a few basic transforms and using the properties allows a large class of signals to 
be easily studied. Examples of modulation, sampling, and others will be given. 


¢ If x(t) = d(t) then X(w) =1 

¢ If x(t) = 1 then X(w) = 27d(w) 

¢ If x(t) is an infinite sequence of delta functions spaced T' apart, 
z(t)=>>~ . 6(¢ — nT), its transform is also an infinite sequence of delta 


functions of weight 27/T spaced 27/T apart, X (w) = 20707. 6 (w — 2ak/T). 
e Other interesting and illustrative examples can be found in [link], [link]. 


Note the Fourier transform takes a function of continuous time into a function of 
continuous frequency, neither function being periodic. If “distribution" or “delta 
functions" are allowed, the Fourier transform of a periodic function will be a infinitely 
long string of delta functions with weights that are the Fourier series coefficients. 


The Laplace Transform 


The Laplace transform can be thought of as a generalization of the Fourier transform in 
order to include a larger class of functions, to allow the use of complex variable theory, to 
solve initial value differential equations, and to give a tool for input-output description of 
linear systems. Its use in system and signal analysis became popular in the 1950's and 
remains as the central tool for much of continuous time system theory. The question of 
convergence becomes still more complicated and depends on complex values of s used in 
the inverse transform which must be in a “region of convergence" (ROC). 


Definition of the Laplace Transform 


The definition of the Laplace transform (LT) of a real valued function defined over all 
positive time tf is 


Equation: 
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and the inverse transform (ILT) is given by the complex contour integral 
Equation: 
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where s = o + jw is a complex variable and the path of integration for the ILT must be in 
the region of the s plane where the Laplace transform integral converges. This definition is 
often called the bilateral Laplace transform to distinguish it from the unilateral transform 
(ULT) which is defined with zero as the lower limit of the forward transform integral 
[link]. Unless stated otherwise, we will be using the bilateral transform. 


Notice that the Laplace transform becomes the Fourier transform on the imaginary axis, 
for s = jw. If the ROC includes the jw axis, the Fourier transform exists but if it does not, 
only the Laplace transform of the function exists. 


There is a considerable literature on the Laplace transform and its use in continuous-time 
system theory. We will develop most of these ideas for the discrete-time system in terms 
of the z-transform later in this chapter and will only briefly consider only the more 
important properties here. 


The unilateral Laplace transform cannot be used if useful parts of the signal exists for 
negative time. It does not reduce to the Fourier transform for signals that exist for negative 
time, but if the negative time part of a signal can be neglected, the unilateral transform will 
converge for a much larger class of function that the bilateral transform will. It also makes 
the solution of linear, constant coefficient differential equations with initial conditions 
much easier. 


Properties of the Laplace Transform 


Many of the properties of the Laplace transform are similar to those for Fourier transform 
[link], [link], however, the basis functions for the Laplace transform are not orthogonal. 
Some of the more important ones are: 


1. Linear: 2{a2 + y} = 2{x} + L{y} 


2. Convolution: If y(t) = h(t)*a(t) = [h(t — 7) x(r) dr 
then #{h(t)*x(t)} = P{h(t)} L{x(t)} 

3. Derivative: {a} =27 {a (t)} 

4. Derivative (ULT): Y { @} = s¥ {x (t)} — x(0) 

5. Integral: Y { f x(t) dt} = +L {x (t)} 

6. Shift: & {x (t — T)} = C(k)e-7s 

7. Modulate: # {a (t) e#"} = X (s — jwo) 


Examples can be found in [link], [link] and are similar to those of the z-transform 
presented later in these notes. Indeed, note the parallals and differences in the Fourier 
series, Fourier transform, and Z-transform. 


Discrete-Time Signals 


Although the discrete-time signal x(n) could be any ordered sequence of numbers, 
they are usually samples of a continuous-time signal. In this case, the real or imaginary 
valued mathematical function x(7) of the integer 7 is not used as an analogy of a 
physical signal, but as some representation of it (such as samples). In some cases, the 
term digital signal is used interchangeably with discrete-time signal, or the label digital 
signal may be use if the function is not real valued but takes values consistent with 
some hardware system. 


Indeed, our very use of the term “discrete-time" indicates the probable origin of the 
signals when, in fact, the independent variable could be length or any other variable or 
simply an ordering index. The term “digital” indicates the signal is probably going to 
be created, processed, or stored using digital hardware. As in the continuous-time case, 
the Fourier transform will again be our primary tool [link], [link], [link]. 


Notation has been an important element in mathematics. In some cases, discrete-time 
signals are best denoted as a sequence of values, in other cases, a vector is created with 
elements which are the sequence values. In still other cases, a polynomial is formed 
with the sequence values as coefficients for a complex variable. The vector formulation 
allows the use of linear algebra and the polynomial formulation allows the use of 
complex variable theory. 


The Discrete Fourier Transform 


The description of signals in terms of their sinusoidal frequency content has proven to 
be as powerful and informative for discrete-time signals as it has for continuous-time 
signals. It is also probably the most powerful computational tool we will use. We now 
develop the basic discrete-time methods starting with the discrete Fourier transform 
(DFT) applied to finite length signals, followed by the discrete-time Fourier transform 
(DTFT) for infinitely long signals, and ending with the z-transform which uses the 
powerful tools of complex variable theory. 


Definition of the DFT 


It is assumed that the signal z(n) to be analyzed is a sequence of N real or complex 
values which are a function of the integer variable n. The DFT of x(n), also called the 
spectrum of x(n), is a length N sequence of complex numbers denoted C'(k) and 
defined by 

Equation: 


N-1 


C(k) = $< a(n) tan 


n= 


using the usual engineering notation: 7 = /—1. The inverse transform (IDFT) which 
retrieves z(n) from C(k) is given by 
Equation: 


N k=0 


C (k) ef ¥"* 


a(n) = 


which is easily verified by substitution into [link]. Indeed, this verification will require 
using the orthogonality of the basis function of the DFT which is 
Equation: 

a N ifn=m 
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k=0 


The exponential basis functions, e irk for k € {0, N — 1}, are the N values of the 
Nth roots of unity (the N zeros of the polynomial (s — 8 ae ). This property is what 
connects the DFT to convolution and allows efficient algorithms for calculation to be 
developed [link]. They are used so often that the following notation is defined by 
Equation: 
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with the subscript being omitted if the sequence length is obvious from context. Using 
this notation, the DFT becomes 


Equation: 
N-1 
C(k) = So a(n) Wyk 
n=0 


One should notice that with the finite summation of the DFT, there is no question of 
convergence or of the ability to interchange the order of summation. No “delta 
functions” are needed and the N transform values can be calculated exactly (within the 


accuracy of the computer or calculator used) from the N signal values with a finite 
number of arithmetic operations. 


Matrix Formulation of the DFT 


There are several advantages to using a matrix formulation of the DFT. This is given by 
writing [link] or [link] in matrix operator form as 


Equation: 
Co Ww? w® W? ... w°® Lo 
Ci Ww? W! W? L1 
C> _ Ww? Ww? w4 “v2 
Cn-1 w? — WN-1)(N-1) ty-1 
or 
Equation: 
C= Fx, 


The orthogonality of the basis function in [link] shows up in this matrix formulation by 
the columns of F’ being orthogonal to each other as are the rows. This means that 
FIF = KI, where k isa scalar constant, and, therefore, F? = kF~—!. This is called 
a unitary operator. 


The definition of the DFT in [link] emphasizes the fact that each of the NV DFT values 
are the sum of NV products. The matrix formulation in [link] has two interpretations. 
Each k-th DFT term is the inner product of two vectors, k-th row of F and x; or, the 
DFT vector, C is a weighted sum of the N columns of F with weights being the 
elements of the signal vector x. A third view of the DFT is the operator view which is 
simply the single matrix equation [link]. 


It is instructive at this point to write a computer program to calculate the DFT of a 
signal. In Matlab [link], there is a pre-programmed function to calculate the DFT, but 
that hides the scalar operations. One should program the transform in the scalar 
interpretive language of Matlab or some other lower level language such as 
FORTRAN, C, BASIC, Pascal, etc. This will illustrate how many multiplications and 


additions and trigonometric evaluations are required and how much memory is needed. 
Do not use a complex data type which also hides arithmetic, but use Euler's relations 
Equation: 

}x 


e/” = cos (x) + 7 sin (2) 


to explicitly calculate the real and imaginary part of C(k). 


If Matlab is available, first program the DFT using only scalar operations. It will 
require two nested loops and will run rather slowly because the execution of loops is 
interpreted. Next, program it using vector inner products to calculate each C(k) which 
will require only one loop and will run faster. Finally, program it using a single matrix 
multiplication requiring no loops and running much faster. Check the memory 
requirements of the three approaches. 


The DFT and IDFT are a completely well-defined, legitimate transform pair with a 
sound theoretical basis that do not need to be derived from or interpreted as an 
approximation to the continuous-time Fourier series or integral. The discrete-time and 
continuous-time transforms and other tools are related and have parallel properties, but 
neither depends on the other. 


The notation used here is consistent with most of the literature and with the standards 
given in [link]. The independent index variable n of the signal x(7) is an integer, but it 
is usually interpreted as time or, occasionally, as distance. The independent index 
variable k of the DFT C(k) is also an integer, but it is generally considered as 
frequency. The DFT is called the spectrum of the signal and the magnitude of the 
complex valued DFT is called the magnitude of that spectrum and the angle or 
argument is called the phase. 


Extensions of x(n) 


Although the finite length signal x(n) is defined only over the interval 

{0 <n < (N — 1)}, the IDFT of C(k) can be evaluated outside this interval to give 
well defined values. Indeed, this process gives the periodic property 4. There are two 
ways of formulating this phenomenon. One is to periodically extend x(n) to —oo and 
+oo and work with this new signal. A second more general way is evaluate all indices 
n and k modulo N. Rather than considering the periodic extension of x(7) on the line 
of integers, the finite length line is formed into a circle or a line around a cylinder so 
that after counting to N — 1, the next number is zero, not a periodic replication of it. 
The periodic extension is easier to visualize initially and is more commonly used for 


the definition of the DFT, but the evaluation of the indices by residue reduction modulo 
N is amore general definition and can be better utilized to develop efficient algorithms 
for calculating the DFT [link]. 


Since the indices are evaluated only over the basic interval, any values could be 
assigned (7) outside that interval. The periodic extension is the choice most 
consistent with the other properties of the transform, however, it could be assigned to 
zero [link]. An interesting possibility is to artificially create a length 2NV sequence by 
appending x(—7) to the end of x(n). This would remove the discontinuities of 
periodic extensions of this new length 2N signal and perhaps give a more accurate 
measure of the frequency content of the signal with no artifacts caused by “end 
effects". Indeed, this modification of the DFT gives what is called the discrete cosine 
transform (DCT) [link]. We will assume the implicit periodic extensions to x(n) with 
no special notation unless this characteristic is important, then we will use the notation 


Convolution 


Convolution is an important operation in signal processing that is in some ways more 
complicated in discrete-time signal processing than in continuous-time signal 
processing and in other ways easier. The basic input-output relation for a discrete-time 
system is given by so-called linear or non-cyclic convolution defined and denoted by 
Equation: 


where x(7) is the perhaps infinitely long input discrete-time signal, h(m) is the 
perhaps infinitely long impulse response of the system, and y(7) is the output. The 
DFT is, however, intimately related to cyclic convolution, not non-cyclic convolution. 
Cyclic convolution is defined and denoted by 

Equation: 


where either all of the indices or independent integer variables are evaluated modulo N 
or all of the signals are periodically extended outside their length N domains. 


This cyclic (sometimes called circular) convolution can be expressed as a matrix 
operation by converting the signal h(n) into a matrix operator as 


Equation: 
ho Ap-1 hyo --: hy 
hy ho At-a 
H—= fro mm ho 
Ries ee; ho 


The cyclic convolution can then be written in matrix notation as 
Equation: 


Y = HX 


where X and Y are column matrices or vectors of the input and output values 
respectively. 


Because non-cyclic convolution is often what you want to do and cyclic convolution is 
what is related to the powerful DFT, we want to develop a way of doing non-cyclic 
convolution by doing cyclic convolution. 


The convolution of a length NV sequence with a length M sequence yields a length 
N + M — 1 output sequence. The calculation of non-cyclic convolution by using 
cyclic convolution requires modifying the signals by appending zeros to them. This 
will be developed later. 


Properties of the DFT 


The properties of the DFT are extremely important in applying it to signal analysis and 
to interpreting it. The main properties are given here using the notation that the DFT of 
a length-N complex sequence x(n) is A{x(n)} = C(k). 


1. Linear Operator: F{z(n) + y(n)} = F{x(n)}+ F{y(n)} 
2. Unitary Operator: F~? = 
3. Periodic Spectrum: C(k) = C(k+N 

4. Periodic Extensions of x(n): z(n) = a(n+N) 


2|- 


De 


WO OND 


13. 
14. 


15. 


Properties of Even and Odd Parts: 2(n) = u(n) + ju(n) and 
C(k) = A(k) + jB(k) 


u v A B |C| 6 
even 0 even 0 even 0 
odd 0 0 odd even m/2 
0 even 0 even even 1/2 
0 odd odd 0 even 0 


. Cyclic Convolution: A{h(n) ox(n)} = F{h(n)}F{ax(n)} 
. Multiplication: F{h(n)z(n)} = F{h(n)}o F{x(n)} 
. Parseval: Yn |@(n)|” = a7 Deco IC (A) 
. Shift: F {a(n —M)} = C(k)e227Mk/N 
. Modulate: F {x (n)e®"™*"/N} = C(k—K) 
. Down Sample or Decimate: ¥ {a (Kn)} = + ee C (k + Lm) where 
N=LK 
. Up Sample or Stretch: If x, (2n) = x(n) for integer n and zero otherwise, then 
Ff{xe;(n)} = C(k), fork =0,1,2,...,2N—1 
N Roots of Unity: (wky” = Lfork = 0,1,2,...,N—1 
Orthogonality: 
Equation: 
= eo dammk/N ,j2mnk/N 9 ifn =m 
—~ 0 ifnA#~m. 
Diagonalization of Convolution: If cyclic convolution is expressed as a matrix 


operation by y = Hx with H given by [link], the DFT operator diagonalizes the 
convolution operator H, or FTHF = Hg where Hg isa diagonal matrix with 
the N values of the DFT of h(n) on the diagonal. This is a matrix statement of 
Property 6. Note the columns of F are the N eigenvectors of H, independent of 
the values of h(n). 


One can show that any “kernel" of a transform that would support cyclic, length-N 
convolution must be the N roots of unity. This says the DFT is the only transform over 
the complex number field that will support convolution. However, if one considers 
various finite fields or rings, an interesting transform, called the Number Theoretic 
Transform, can be defined and used because the roots of unity are simply two raised to 
a powers which is a simple word shift for certain binary number representations [link], 
[link]. 


Examples of the DFT 


It is very important to develop insight and intuition into the DFT or spectral 
characteristics of various standard signals. A few DFT's of standard signals together 
with the above properties will give a fairly large set of results. They will also aid in 
quickly obtaining the DFT of new signals. The discrete-time impulse 6(7) is defined 
by 

Equation: 


5(n) _ 1 va 
0 otherwise 


The discrete-time pulse [yy (7) is defined by 


Equation: 
1 when n=0,1,---,M-1 
0 otherwise 


Om(n) = { 


Several examples are: 


¢ DFT{6(n)} = 1, The DFT of an impulse is a constant. 

¢ DFT{1} = N6(k), The DFT of a constant is an impulse. 

© DFT {e®*Kn/N\ — N6(k— K) 

¢ DFT {cos (2nMn/N) = 1[6(k—- M)+6(k+M)] 
k) 


° DFT {ny (n)} = a - 


These examples together with the properties can generate a still larger set of interesting 
and enlightening examples. Matlab can be used to experiment with these results and to 
gain insight and intuition. 


The Discrete-Time Fourier Transform 


In addition to finite length signals, there are many practical problems where we must 
be able to analyze and process essentially infinitely long sequences. For continuous- 
time signals, the Fourier series is used for finite length signals and the Fourier 
transform or integral is used for infinitely long signals. For discrete-time signals, we 
have the DFT for finite length signals and we now present the discrete-time Fourier 
transform (DTFT) for infinitely long signals or signals that are longer than we want to 
specify [link]. The DTFT can be developed as an extension of the DFT as N goes to 
infinity or the DTFT can be independently defined and then the DFT shown to be a 
special case of it. We will do the latter. 


Definition of the DTFT 


The DTFT of a possibly infinitely long real (or complex) valued sequence f(n) is 
defined to be 
Equation: 


and its inverse denoted IDTFT is given by 
Equation: 


f(n) = = 7 F (w) e%" dw. 


Verification by substitution is more difficult than for the DFT. Here convergence and 
the interchange of order of the sum and integral are serious questions and have been the 
topics of research over many years. Discussions of the Fourier transform and series for 
engineering applications can be found in [link], [link]. It is necessary to allow 
distributions or delta functions to be used to gain the full benefit of the Fourier 
transform. 


Note that the definition of the DIFT and IDTFT are the same as the definition of the 
IFS and FS respectively. Since the DTFT is a continuous periodic function of w, its 
Fourier series is a discrete set of values which turn out to be the original signal. This 
duality can be helpful in developing properties and gaining insight into various 
problems. The conditions on a function to determine if it can be expanded in a FS are 


exactly the conditions on a desired frequency response or spectrum that will determine 
if a signal exists to realize or approximate it. 


Properties 


The properties of the DTFT are similar to those for the DFT and are important in the 
analysis and interpretation of long signals. The main properties are given here using the 
notation that the DTFT of a complex sequence x(n) is F{a(n)} = X(w). 


1. Linear Operator: F{x+y} = ¥{«x}+ F{y} 

2. Periodic Spectrum: X(w) = X(w + 27) 

3. Properties of Even and Odd Parts: x(n) = u(n) + ju(n) and 
X(w) = A(w) + jBW) 


U Vv A B |X| 6 
even 0 even 0 even 0 
odd 0 0 odd even 0 
0 even 0 even even m/2 
0 odd odd 0 even m/2 


4. Convolution: If non-cyclic or linear convolution is defined by: 
y(n) = h(n)*a(n) = Yon h(n — m)a(m) = die. h(k)a (n — k) 
then F{h(n)*x(n)} = F{h(n)}.F{x(n)} 
5. Multiplication: If cyclic convolution is defined by: 
Y(w) = H(w)oX(w) = fo Hw—2)X (da 
F {h(n)x(n)} = a¢F {h(n)} oF {x (n)} 
6. Parseval: 7°24 |x (n)/? = 2 f™ |X (w)|?dw 
7. Shift: F {a(n —M)} = X (w)e-*™ 
8. Modulate: F {x (n)e"} = X(w—wy) 
9. Sample: ¥ {x(Kn)} = 4 Saw X (w+ Lm) where N = LK 


10. Stretch: F{x,(n)} = X (w), for -—Ka <w < Knwherex,(Kn) = x(n) 
for integer n and zero otherwise. — 
11. Orthogonality: \>°° eT" EI" = 27d (wy — we) 


Evaluation of the DTFT by the DFT 


If the DTFT of a finite sequence is taken, the result is a continuous function of w. If the 
DFT of the same sequence is taken, the results are N evenly spaced samples of the 
DTFT. In other words, the DTFT of a finite signal can be evaluated at N points with 
the DFT. 

Equation: 


X(w) = DTFT {x(n)} = S- a (n)e 


and because of the finite length 
Equation: 


N-1 


X(w) = S> a(n) eI", 


n= 


If we evaluate w at NV equally space points, this becomes 
Equation: 


which is the DFT of x(n). By adding zeros to the end of x(n) and taking a longer 
DFT, any density of points can be evaluated. This is useful in interpolation and in 
plotting the spectrum of a finite length signal. This is discussed further in Sampling, 
Up-Sampling, Down-Sampling, and Multi-Rate Processing. 


There is an interesting variation of the Parseval's theorem for the DTFT of a finite 
length-N signal. If x(n) 4 0 for0 > n > N — 1, andif L > N, then 
Equation: 


2» 1 2 : 2 
z(n)|° = — X (27k/L)\" = — X (w)|" dw. 
Sey? = FY x eak/E)? = = fx) 


n= 


The second term in [link] says the Riemann sum is equal to its limit in this case. 


Examples of DTFT 


As was true for the DFT, insight and intuition is developed by understanding the 
properties and a few examples of the DTFT. Several examples are given below and 
more can be found in the literature [link], [link], [link]. Remember that while in the 
case of the DFT signals were defined on the region {0 < n < (N — 1)} and values 
outside that region were periodic extensions, here the signals are defined over all 
integers and are not periodic unless explicitly stated. The spectrum is periodic with 
period 27. 


¢ DTFT{6(n)} = 1 for all frequencies. 
e Equation: 


DTFT{1} = 27d(w) 
e Equation: 
DTFT {e%"} = 276 (w — wo) 


e Equation: 


DTFT {cos (won)} = m[d(w— wo) + 6 (w+ wo)] 
e Equation: 


sin (wMk/2) 


DTFT {1m (n)} = sin (wk/2) 


The Z-Transform 


The z-transform is an extension of the DTFT in a way that is analogous to the Laplace 
transform for continuous-time signals being an extension of the Fourier transform. It 
allows the use of complex variable theory and is particularly useful in analyzing and 
describing systems. The question of convergence becomes still more complicated and 


depends on values of z used in the inverse transform which must be in the “region of 
convergence" (ROC). 


Definition of the Z-Transform 


The z-transform (ZT) is defined as a polynomial in the complex variable z with the 
discrete-time signal values as its coefficients [link], [link], [link]. It is given by 
Equation: 


CO 
F(z) = >> f(n)z™ 
n=—oo 
and the inverse transform (IZT) is 
Equation: 
1 n—-1 
f(n) = — F(z) z” “dz. 
27j J ROC 


The inverse transform can be derived by using the residue theorem [link], [link] from 
complex variable theory to find f(0) from z~! F(z), f(1) from F(z), f(2) from 
zF(z), and in general, f(n) from z”~!F (z). Verification by substitution is more 
difficult than for the DFT or DTFT. Here convergence and the interchange of order of 
the sum and integral is a serious question that involves values of the complex variable 
z. The complex contour integral in [link] must be taken in the ROC of the z plane. 


A unilateral z-transform is sometimes needed where the definition [link] uses a lower 
limit on the transform summation of zero. This allow the transformation to converge 
for some functions where the regular bilateral transform does not, it provides a 
straightforward way to solve initial condition difference equation problems, and it 
simplifies the question of finding the ROC. The bilateral z-transform is used more for 
signal analysis and the unilateral transform is used more for system description and 
analysis. Unless stated otherwise, we will be using the bilateral z-transform. 


Properties 


The properties of the ZT are similar to those for the DTFT and DFT and are important 
in the analysis and interpretation of long signals and in the analysis and description of 


discrete-time systems. The main properties are given here using the notation that the 
ZT of a complex sequence x(n) is &{x(n)} = X(z). 


1 
2 
3: 
4 


uo 


fp) 


N 


jee) 


ice) 


. Modulate: #{z (n)a”} 
10. 
ine 


. Linear Operator: @{x+y} = #{x}+ #{y} 
. Relationship of ZT to DTFT: & {x}|,_.4 = DZ FT {x} 
Periodic Spectrum: X (e”) = X (e+?) 
. Properties of Even and Odd Parts: x(n) = u(n) + ju(n) and 
X(e#) = A(e*) + jB(e*) 
Equation: 
U Vv A B 
even QO even 0 
odd 0 0 odd 
0 even O- even 
0 odd_ odd 0 
. Convolution: If discrete non-cyclic convolution is defined by 


y(n) = h(n)*a(n) = Yim=—oo h(n —m)x(m) = Dip h (k)a (n — k) 
then #{h(n)*ax(n)} = #{h(n)}#{a(n)} 


. Shift: # {2 (n+ M)} = 2”@X (2) 
. Shift (unilateral): 

File(n+m)} = 2"X (z) —2"2 (0) —2" 12 (1) —--- — za (m—1) 
. Shift (unilateral): 


#¥{e(n—m)} = 2-™X (z) — 2-™ 12 (-1) —--- — 2 (—m) 
= X(z/a) 
Time mult.: ¥ {nz (n)} = (—z)™ ie 
Evaluation: The ZT can be evaluated on the unit circle in the z-plane by taking the 
DTFT of x(n) and if the signal is finite in length, this can be evaluated at sample 


points by the DFT. 


Examples of the Z-Transform 


A few examples together with the above properties will enable one to solve and 
understand a wide variety of problems. These use the unit step function to remove the 
negative time part of the signal. This function is defined as 

Equation: 


(n) 1 ifn>0 
u(n) = 
0 ifn<0 


and several bilateral z-transforms are given by 


© £{d(n)} = 1 forall z. 
\ = tor |4| >, 


(n)} = 
© #{u(n)a"} = — for |z| > |al. 


Z—4 


Notice that these are similar to but not the same as a term of a partial fraction 
expansion. 


Inversion of the Z-Transform 


The z-transform can be inverted in three ways. The first two have similar procedures 
with Laplace transformations and the third has no counter part. 


e The z-transform can be inverted by the defined contour integral in the ROC of the 
complex z plane. This integral can be evaluated using the residue theorem [link], 
[link]. 

e The z-transform can be inverted by expanding iF (z) ina partial fraction 
expansion followed by use of tables for the first or second order terms. 


AS , f(n) can be 


e The third method is not analytical but numerical. If F(z) = 
obtained as the coefficients of long division. 
For example 
Equation: 


=14tazr t@rs? i+: 


z—a 


which is u (7) a” as used in the examples above. 


We must understand the role of the ROC in the convergence and inversion of the z- 
transform. We must also see the difference between the one-sided and two-sided 
transform. 


Solution of Difference Equations using the Z-Transform 


The z-transform can be used to convert a difference equation into an algebraic equation 
in the same manner that the Laplace converts a differential equation in to an algebraic 
equation. The one-sided transform is particularly well suited for solving initial 
condition problems. The two unilateral shift properties explicitly use the initial values 
of the unknown variable. 


A difference equation DE contains the unknown function x(n) and shifted versions of 
it such as z(n — 1) or x(n + 3). The solution of the equation is the determination of 
z(t). A linear DE has only simple linear combinations of x(n) and its shifts. An 
example of a linear second order DE is 

Equation: 


ax(n)+ba(n—1)+ca(n—2) = f(n) 


A time invariant or index invariant DE requires the coefficients not be a function of n 
and the linearity requires that they not be a function of x(n). Therefore, the 
coefficients are constants. 


This equation can be analyzed using classical methods completely analogous to those 
used with differential equations. A solution of the form x(n) = KX” is substituted 
into the homogeneous difference equation resulting in a second order characteristic 
equation whose two roots give a solution of the form a, (n) = KyAT + KoAf. A 
particular solution of a form determined by f(n) is found by the method of 
undetermined coefficients, convolution or some other means. The total solution is the 
particular solution plus the solution of the homogeneous equation and the three 
unknown constants K; are determined from three initial conditions on x(n). 


It is possible to solve this difference equation using z-transforms in a similar way to the 
solving of a differential equation by use of the Laplace transform. The z-transform 
converts the difference equation into an algebraic equation. Taking the ZT of both sides 
of the DE gives 

Equation: 


aX (z)+b[z7X (z)+2(-1)] +e [z 7X (2) +22 (-1) +. z(-2)] = Y(z) 


solving for X(z) gives 
Equation: 


27 [Y (z) — ba (-1) — x(-2)] — zex(-1) 
az?+bz+c 


and inversion of this transform gives the solution x(n). Notice that two initial values 
were required to give a unique solution just as the classical method needed two values. 


These are very general methods. To solve an nth order DE requires only factoring ann 
th order polynomial and performing a partial fraction expansion, jobs that computers 
are well suited to. There are problems that crop up if the denominator polynomial has 
repeated roots or if the transform of y(n) has a root that is the same as the 
homogeneous equation, but those can be handled with slight modifications giving 
solutions with terms of the from nA” just as similar problems gave solutions for 
differential equations of the form t e*’. 


The original DE could be rewritten in a different form by shifting the index to give 
Equation: 


axz(n+2)+ba(n+1)+ca(n) = f(n+2) 


which can be solved using the second form of the unilateral z-transform shift property. 


Region of Convergence for the Z-Transform 


Since the inversion integral must be taken in the ROC of the transform, it is necessary 
to understand how this region is determined and what it means even if the inversion is 
done by partial fraction expansion or long division. Since all signals created by linear 
constant coefficient difference equations are sums of geometric sequences (or samples 
of exponentials), an analysis of these cases will cover most practical situations. 
Consider a geometric sequence starting at zero. 

Equation: 


with a z-transform 
Equation: 


F(z) = 14t@e7) bee tae te eae. 


Multiplying by a z~! gives 
Equation: 


az (FQ) = ad? fas? +a 7 tae 4 ta" 


and subtracting from [link] gives 
Equation: 


(l—az")F(z) = 1- i mes 


Solving for F(z) results in 


Equation: 
1 — gMti,-M-1 z—a(2)™ 
F(z) = —__—_ = 
1—az1 2=G 


The limit of this sum as MZ — ov is 
Equation: 


for |z| > |a|. This not only establishes the z-transform of f(n) but gives the region in 
the z plane where the sum converges. 


If a similar set of operations is performed on the sequence that exists for negative n 
Equation: 

a” n<0 

0 n>0O 


the result is 
Equation: 


for |z| < |a|. Here we have exactly the same z-transform for a different sequence f(n) 
but with a different ROC. The pole in F(z) divides the z-plane into two regions that 
give two different f(n). This is a general result that can be applied to a general rational 
F(z) with several poles and zeros. The z-plane will be divided into concentric annular 
regions separated by the poles. The contour integral is evaluated in one of these regions 
and the poles inside the contour give the part of the solution existing for negative n 
with the poles outside the contour giving the part of the solution existing for positive n. 


Notice that any finite length signal has a z-transform that converges for all z. The ROC 
is the entire z-plane except perhaps zero and/or infinity. 


Relation of the Z-Transform to the DTFT and the DFT 


The FS coefficients are weights on the delta functions in a FT of the periodically 
extended signal. The FT is the LT evaluated on the imaginary axis:s = jw. 


The DFT values are samples of the DTFT of a finite length signal. The DTFT is the z- 
transform evaluated on the unit circle in the z plane. 


Equation: 
Bo). = S> eile ! =2 740m} 

Equation: 

Ble) = S° a(nje I" = GIFT {x(n)} 
and if x(n) is of length NV 
Equation: 

. an Wee . an 
F (ei**) = a(njeIv" — GFF {x(n)} 
n=0 


It is important to be able to relate the time-domain signal x(n), its spectrum X(w), and 
its z-transform represented by the pole-zero locations on the z plane. 


Relationships Among Fourier Transforms 


The DFT takes a periodic discrete-time signal into a periodic discrete-frequency 
representation. 


The DTFT takes a discrete-time signal into a periodic continuous-frequency 
representation. 


The FS takes a periodic continuous-time signal into a discrete-frequency 
representation. 


The FT takes a continuous-time signal into a continuous-frequency representation. 


The LT takes a continuous-time signal into a function of a continuous complex 
variable. 


The ZT takes a discrete-time signal into a function of a continuous complex variable. 


Wavelet-Based Signal Analysis 


There are wavelet systems and transforms analogous to the DFT, Fourier series, 
discrete-time Fourier transform, and the Fourier integral. We will start with the discrete 
wavelet transform (DWT) which is analogous to the Fourier series and probably should 
be called the wavelet series [link]. Wavelet analysis can be a form of time-frequency 
analysis which locates energy or events in time and frequency (or scale) 
simultaneously. It is somewhat similar to what is called a short-time Fourier transform 
or a Gabor transform or a windowed Fourier transform. 


The history of wavelets and wavelet based signal processing is fairly recent. Its roots in 
signal expansion go back to early geophysical and image processing methods and in 
DSP to filter bank theory and subband coding. The current high interest probably 
started in the late 1980's with the work of Mallat, Daubechies, and others. Since then, 
the amount of research, publication, and application has exploded. Two excellent 
descriptions of the history of wavelet research and development are by Hubbard [link] 
and by Daubechies [link] and a projection into the future by Sweldens [link] and 
Burrus [link]. 


The Basic Wavelet Theory 


The ideas and foundations of the basic dyadic, multiresolution wavelet systems are 
now pretty well developed, understood, and available [link], [link], [link], [link]. The 
first basic requirement is that a set of expansion functions (usually a basis) are 
generated from a single “mother” function by translation and scaling. For the discrete 
wavelet expansion system, this is 


Equation: 


yj,n (t) = y (2’t — k) 


where 7, & are integer indices for the series expansion of the form 
Equation: 


f(t) = do cin p5a (2). 
j,k 


The coefficients c;, are called the discrete wavelet transform of the signal f(t). This 


use of translation and scale to create an expansion system is the foundation of all so- 
called first generation wavelets [Link]. 


The system is somewhat similar to the Fourier series described in Equation 51 from 
Least Squared Error Designed of FIR Filters with frequencies being related by powers 
of two rather than an integer multiple and the translation by k giving only the two 
results of cosine and sine for the Fourier series. 


The second almost universal requirement is that the wavelet system generates a 
multiresolution analysis (MRA). This means that a low resolution function (low scale 7 
) can be expanded in terms of the same function at a higher resolution (higher 7). This 
is stated by requiring that the generator of a MRA wavelet system, called a scaling 
function y(t), satisfies 

Equation: 


9 (t) = h(n) 9 (2t —n). 


This equation, called the refinement equation or the MRA equation or basic 
recursion equation, is similar to a differential equation in that its solution is what 
defines the basic scaling function and wavelet [link], [link]. 


The current state of the art is that most of the necessary and sufficient conditions on the 
coefficients h(n) are known for the existence, uniqueness, orthogonality, and other 
properties of y(t). Some of the theory parallels Fourier theory and some does not. 


A third important feature of a MRA wavelet system is a discrete wavelet transform 
(DWT) can be calculated by a digital filter bank using what is now called Mallat's 
algorithm. Indeed, this connection with digital signal processing (DSP) has been a rich 
source of ideas and methods. With this filter bank, one can calculate the DWT of a 


length-N digital signal with order N operations. This means the number of 
multiplications and additions grows only linearly with the length of the signal. This 
compares with N log (NV) for an FFT and N? for most methods and worse than that 
for some others. 


These basic ideas came from the work of Meyer, Daubechies, Mallat, and others but 
for a time looked like a solution looking for a problem. Then a second phase of 
research showed there are many problems to which the wavelet is an excellent solution. 
In particular, the results of Donoho, Johnstone, Coifman, Beylkin, and others opened 
another set of doors. 


Generalization of the Basic Wavelet System 


After (in some cases during) much of the development of the above basic ideas, a 
number of generalizations [link] were made. They are listed below: 


1. A larger integer scale factor than MZ = 2 can be used to give a more general M- 
band refinement equation [link] 
Equation: 


p(t) = (mn) p (Mt — n) 


than the “dyadic" or octave based Equation 4 from Rational Function 
Approximation. This also gives more than two channels in the accompanying 
filter bank. It allows a uniform frequency resolution rather than the resulting 
logarithmic one for M = 2. 

. The wavelet system called a wavelet packet is generated by “iterating” the 
wavelet branches of the filter bank to give a finer resolution to the wavelet 
decomposition. This was suggested by Coifman and it too allows a mixture of 
uniform and logarithmic frequency resolution. It also allows a relatively simple 
adaptive system to be developed which has an automatically adjustable frequency 
resolution based on the properties of the signal. 

. The usual requirement of translation orthogonality of the scaling function and 
wavelets can be relaxed to give what is called a biorthogonal system{[link]. If the 
expansion basis is not orthogonal, a dual basis can be created that will allow the 
usual expansion and coefficient calculations to be made. The main disadvantage is 
the loss of a Parseval's theorem which maintains energy partitioning. 
Nevertheless, the greater flexibility of the biorthogonal system allows superior 
performance in many compression and denoising applications. 

4. The basic refinement Equation 4 from Rational Function Approximation gives the 

scaling function in terms of a compressed version of itself (self-similar). If we 


NO 


ise) 


allow two (or more) scaling functions, each being a weighted sum of a compress 
version of both, a more general set of basis functions results. This can be viewed 
as a vector of scaling functions with the coefficients being a matrix now. Once 
again, this generalization allows more flexibility in the characteristics of the 
individual scaling functions and their related multi-wavelets. These are called 
multi-wavelet systems and are still being developed. 

5. One of the very few disadvantages of the discrete wavelet transform is the fact it 
is not shift invariant. In other words, if you shift a signal in time, its wavelet 
transform not only shifts, it changes character! For many applications in denoising 
and compression, this is not desirable although it may be tolerable. The DWT can 
be made shift-invariant by calculating the DWT of a signal for all possible shifts 
and adding (or averaging) the results. That turns out to be equivalent to removing 
all of the down-samplers in the associated filter bank (an undecimated filter 
bank), which is also equivalent to building an overdetermined or redundant 
DWT from a traditional wavelet basis. This overcomplete system is similar to a 
“tight frame" and maintains most of the features of an orthogonal basis yet is shift 
invariant. It does, however, require N log (V) operations. 

6. Wavelet systems are easily modified to being an adaptive system where the basis 
adjusts itself to the properties of the signal or the signal class. This is often done 
by starting with a large collection or library of expansion systems and bases. A 
subset is adaptively selected based on the efficiency of the representation using a 
process sometimes called pursuit. In other words, a set is chosen that will result 
in the smallest number of significant expansion coefficients. Clearly, this is signal 
dependent, which is both its strength and its limitation. It is nonlinear. 

7. One of the most powerful structures yet suggested for using wavelets for signal 
processing is to first take the DWT, then do a point-wise linear or nonlinear 
processing of the DWT, finally followed by an inverse DWT. Simply setting some 
of the wavelet domain expansion terms to zero results in linear wavelet domain 
filtering, similar to what would happen if the same were done with Fourier 
transforms. Donoho [link], [link] and others have shown by using some form of 
nonlinear thresholding of the DWT, one can achieve near optimal denoising or 
compression of a signal. The concentrating or localizing character of the DWT 
allows this nonlinear thresholding to be very effective. 


The present state of activity in wavelet research and application shows great promise 
based on the above generalizations and extensions of the basic theory and structure 
[link]. We now have conferences, workshops, articles, newsletters, books, and email 
groups that are moving the state of the art forward. More details, examples, and 
software are given in [link], [link], [link]. 


Discrete-Time Systems 


In the context of discussing signal processing, the most general definition of a 
system is similar to that of a function. A system is a device, formula, rule, or 
some process that assigns an output signal from some given class to each possible 
input signal chosen from some allowed class. From this definition one can pose 
three interesting and practical problems. 


1. Analysis: If the input signal and the system are given, find the output signal. 
2. Control: If the system and the output signal are given, find the input signal. 
3. Synthesis: If the input signal and output signal are given, find the system. 


The definition of input and output signal can be quite diverse. They could be 
scalars, vectors, functions, functionals, or other objects. 


All three of these problems are important, but analysis is probably the most basic 
and its study usually precedes that of the other two. Analysis usually results in a 
unique solution. Control is often unique but there are some problems where 
several inputs would give the same output. Synthesis is seldom unique. There are 
usually many possible systems that will give the same output for a given input. 


In order to develop tools for analysis, control, and design of discrete-time 
systems, specific definitions, restrictions, and classifications must be made. It is 
the explicit statement of what a system is, not what it isn't, that allows a 
descriptive theory and design methods to be developed. 


Classifications 


The basic classifications of signal processing systems are defined and listed here. 
We will restrict ourselves to discrete-time systems that have ordered sequences of 
real or complex numbers as inputs and outputs and will denote the input sequence 
by x(n) and the output sequence by y(m) and show the process of the system by 
x(n) — y(n). Although the independent variable n could represent any physical 
variable, our most common usages causes us to generically call it time but the 
results obtained certainly are not restricted to this interpretation. 


1. Linear, A system is classified as linear if two conditions are true. 


o If z(n) > y(n) thenaaz(n) — ay(n) for all a. This property is called 
homogeneity or scaling. 


o If x1 (n) > yi (n) and x2 (n) > ye (n), then 
(x1 (nm) + £2 (n)) > (y1 (n) + ye (n)) for all x1 and x9. This 
property is called superposition or additivity. 


If a system does not satisfy both of these conditions for all inputs, it is 
classified as nonlinear. For most practical systems, one of these conditions 
implies the other. Note that a linear system must give a zero output for a zero 
input. 

2. Time Invariant , also called index invariant or shift invariant. A system is 
classified as time invariant if z(n + k) — y(n + k) for any integer k. This 
states that the system responds the same way regardless of when the input is 
applied. In most cases, the system itself is not a function of time. 

3. Stable . A system is called bounded-input bounded-output stable if for all 
bounded inputs, the corresponding outputs are bounded. This means that the 
output must remain bounded even for inputs artificially constructed to 
maximize a particular system's output. 

4. Causal . A system is classified as causal if the output of a system does not 
precede the input. For linear systems this means that the impulse response of 
a system is zero for time before the input. This concept implies the 
interpretation of n as time even though it may not be. A system is semi- 
causal if after a finite shift in time, the impulse response is zero for negative 
time. If the impulse response is nonzero for m —> —oo, the system is 
absolutely non-causal. Delays are simple to realize in discrete-time systems 
and semi-causal systems can often be made realizable if a time delay can be 
tolerated. 

5. Real-Time . A discrete-time system can operate in “real-time” if an output 
value in the output sequence can be calculated by the system before the next 
input arrives. If this is not possible, the input and output must be stored in 
blocks and the system operates in “batch" mode. In batch mode, each output 
value can depend on all of the input values and the concept of causality does 
not apply. 


These definitions will allow a powerful class of analysis and design methods to be 
developed and we start with convolution. 


Convolution 


The most basic and powerful operation for linear discrete-time system analysis, 
control, and design is discrete-time convolution. We first define the discrete-time 
unit impulse, also known as the Kronecker delta function, as 


Equation: 


s(n) = 1) for n= 0 


0 otherwise. 


If a system is linear and time-invariant, and 6(n) — h(n), the output y(n) can be 
calculated from its input x(n) by the operation called convolution denoted and 
defined by 

Equation: 


It is informative to methodically develop this equation from the basic properties 
of a linear system. 


Derivation of the Convolution Sum 


We first define a complete set of orthogonal basis functions by 6(n — m) for 

m = 0,1,2,---,00. The input x(n) is broken down into a set of inputs by taking 
an inner product of the input with each of the basis functions. This produces a set 
of input components, each of which is a single impulse weighted by a single value 
of the input sequence (a(n), d(n — m)) = x(m)d(n — m). Using the time 
invariant property of the system, d(n — m) — h(n — m) and using the scaling 
property of a linear system, this gives an output of 

z(m)d(n — M) + «(m)h(n — m). We now calculate the output due to x(n) by 
adding outputs due to each of the resolved inputs using the superposition property 
of linear systems. This is illustrated by the following diagram: 


Equation: 
z(n)d(n) = 2(0)d(n) — 2(0)h(n) 
a(n)d(n—-1) = a(l)d(n-1) 7 a(1)h(n—-1) 

a(n) =< a(n)d(n—2) = x(2)d(n—2) — 2(2)h(n - 2) = y(n) 
sae —m) = a(m)d(n—m) > Se —m) 


or 
Equation: 


CO 


y(n) = 7 2(m)h(n—m) 


m=— Co 


and changing variables gives 
Equation: 


jos Ss eee 


m=— CoO 


If the system is linear but time varying, we denote the response to an impulse at 
n =m by 6(n — m) > h(n,m). In other words, each impulse response may be 
different depending on when the impulse is applied. From the development 
above, it is easy to see where the time-invariant property was used and to derive a 
convolution equation for a time-varying system as 

Equation: 


y(n) =h(n,m)*e(n) = > h(nym)z(m). 


m=—cCo 


Unfortunately, relaxing the linear constraint destroys the basic structure of the 
convolution sum and does not result in anything of this form that is useful. 


By a change of variables, one can easily show that the convolution sum can also 
be written 
Equation: 


y(n) = h(n)*z (n) = 3 h(m)z(n—m). 


m=— CO 


If the system is causal, h(n) = 0 for n < 0 and the upper limit on the summation 
in Equation 2 from Discrete Time Signals becomes m = n. If the input signal is 


causal, the lower limit on the summation becomes zero. The form of the 
convolution sum for a linear, time-invariant, causal discrete-time system with a 
causal input is 

Equation: 


y (n) = h(n)*2 (n) = J> h(n — m)z (m) 
m=0 


or, showing the operations commute 
Equation: 


y(n) =h(n)*x(n) = S° h(m)a(n—m). 
m=0 


Convolution is used analytically to analyze linear systems and it can also be used 
to calculate the output of a system by only knowing its impulse response. This is 
a very powerful tool because it does not require any detailed knowledge of the 
system itself. It only uses one experimentally obtainable response. However, this 
summation cannot only be used to analyze or calculate the response of a given 
system, it can be an implementation of the system. This summation can be 
implemented in hardware or programmed on a computer and become the signal 
processor. 


The Matrix Formulation of Convolution 


Some of the properties and characteristics of convolution and of the systems it 
represents can be better described by a matrix formulation than by the summation 
notation. The first LD values of the discrete-time convolution defined above can be 
written as a matrix operator on a vector of inputs to give a vector of the output 
values. 

Equation: 


U0 ho 0 0 --- 0 LO 


V1 hy ho 0 1 
Y2 — ho hy ho r2 
YL-1 hr-1 vee ho| | ez-1 


If the input sequence z is of length N and the operator signal h is of length M/, 
the output is of length LZ = N + M — 1. This is shown for NV = 4 and M = 3 by 
the rectangular matrix operation 


Equation: 
U0 ho 0 0 0 
Yi hi ho O Of {20 
yo) ho hi ho O L1 
y3| | 0 he hi ho] |ae 
Y4 0 O he hy} [x3 
U5 0 0 0 he 


It is clear that if the system is causal (h(n) = 0 for n < 0), the H matrix is lower 
triangular. It is also easy to see that the system being time-invariant is equivalent 
to the matrix being Toeplitz [link]. This formulation makes it obvious that if a 
certain output were desired from a length 4 input, only 4 of the 6 values could be 
specified and the other 2 would be controlled by them. 


Although the formulation of constructing the matrix from the impulse response of 
the system and having it operate on the input vector seems most natural, the 
matrix could have been formulated from the input and the vector would have been 
the impulse response. Indeed, this might the appropriate formulation if one were 
specifying the input and output and designing the system. 


The basic convolution defined in [link], derived in [link], and given in matrix 
form in [link] relates the input to the output for linear systems. This is the form of 
convolution that is related to multiplication of the DTFT and z-transform of 
signals. However, it is cyclic convolution that is fundamentally related to the DFT 
and that will be efficiently calculated by the fast Fourier transform (FFT) 


developed in Part III of these notes. Matrix formulation of length-L cyclic 
convolution is given by 


Equation: 
Yo ho hea he» -:: hy LO 
Y1 hy ho hy-1 ho Ly 
yo | — | ho hy ho hs x9 
YL-1 hy-1 vee ho} | 1-1 


This matrix description makes it clear that the matrix operator is always square 
and the three signals, x(n), h(n), and y(n), are necessarily of the same length. 


There are several useful conclusions that can be drawn from linear algebra [link]. 
The eigenvalues of the non-cyclic are all the same since the eigenvalues of a 
lower triangular matrix are simply the values on the diagonal. 


Although it is less obvious, the eigenvalues of the cyclic convolution matrix are 
the N values of the DFT of h(n) and the eigenvectors are the basis functions of 
the DFT which are the column vectors of the DFT matrix. The eigenvectors are 
completely controlled by the structure of H being a cyclic convolution matrix and 
are not at all a function of the values of h(n). The DFT matrix equation from 
[link] is given by 

Equation: 


X=Fx and Y=Fy 
where X is the length-N vector of the DFT values, H is the matrix operator for 


the DFT, and x is the length-N vector of the signal z(n) values. The same is true 
for the comparable terms in y. 


The matrix form of the length-N cyclic convolution in [link] is written 
Equation: 


y = Hx 


Taking the DFT both sides and using the IDFT on z gives 


Equation: 


Fy = Y = FHx = FHF 'X 


If we define the diagonal matrix Hg as an L by L matrix with the values of the 
DFT of h(n) on its diagonal, the convolution property of the DFT becomes 
Equation: 


Y = HgX 


This implies 
Equation: 


Hy =FHF! and H=F—'!HgF 


which is the basis of the earlier statement that the eigenvalues of the cyclic 
convolution matrix are the values of the DFT of h(n) and the eigenvectors are the 
orthogonal columns of F’. The DFT matrix diagonalizes the cyclic convolution 
matrix. This is probably the most concise statement of the relation of the DFT to 
convolution and to linear systems. 


An important practical question is how one calculates the non-cyclic convolution 
needed by system analysis using the cyclic convolution of the DFT. The answer is 
easy to see using the matrix description of H. The length of the output of non- 
cyclic convolution is N + M — 1. If N — 1 zeros are appended to the end of 
h(n) and M — 1 zeros are appended to the end of x(n), the cyclic convolution of 
these two augmented signals will produce exactly the same NV + M — 1 values as 
non-cyclic convolution would. This is illustrated for the example considered 
before. 


Equation: 
Yo ho 0 0 O ho hy] 120 
Y1 hy ho O O O hol lay 
yo} [ho hi ho O O OO} Jae 
y3| |0 ho hy ho O OO] Ia 
YA 0 O hog hy ho O 0 


Just enough zeros were appended so that the nonzero terms in the upper right- 
hand corner of H are multiplied by the zeros in the lower part of x and, therefore, 
do not contribute to y. This does require convolving longer signals but the output 
is exactly what we want and we calculated it with the DFT-compatible cyclic 
convolution. Note that more zeros could have been appended to A and = and the 
first NV + M — 1 terms of the output would have been the same only more 
calculations would have been necessary. This is sometimes done in order to use 
forms of the FFT that require that the length be a power of two. 


If fewer zeros or none had been appended to hf and zx, the nonzero terms in the 
upper right-hand corner of H, which are the “tail" of h(n), would have added the 
values that would have been at the end of the non-cyclic output of y(n) to the 
values at the beginning. This is a natural part of cyclic convolution but is 
destructive if non-cyclic convolution is desired and is called aliasing or folding 
for obvious reasons. Aliasing is a phenomenon that occurs in several arenas of 
DSP and the matrix formulation makes it easy to understand. 


The Z-Transform Transfer Function 


Although the time-domain convolution is the most basic relationship of the input 
to the output for linear systems, the z-transform is a close second in importance. It 
gives different insight and a different set of tools for analysis and design of linear 
time-invariant discrete-time systems. 


If our system in linear and time-invariant, we have seen that its output is given by 
convolution. 
Equation: 


= eee 


Mm=— Co 


Assuming that h(7) is such that the summation converges properly, we can 
calculate the output to an input that we already know has a special relation with 
discrete-time transforms. Let x (n) = 2” which gives 

Equation: 


y(n) = \- h(n—m)z™ 


Mm>=— Co 


With the change of variables of k = n — m, we have 


Equation: 
y(n)= So h(k)2e”* =] SO h(k)e*} 2" 
k=—0o k=—0o 
or 
Equation: 


y(n) =H (z)2" 


We have the remarkable result that for an input of x (n) = z”, we get an output 
of exactly the same form but multiplied by a constant that depends on z and this 
constant is the z-transform of the impulse response of the system. In other words, 
if the system is thought of as a matrix or operator, z” is analogous to an 
eigenvector of the system and H(z) is analogous to the corresponding eigenvalue. 


We also know from the properties of the z-transform that convolution in the n 
domain corresponds to multiplication in the z domain. This means that the z- 
transforms of x(n) and y(n) are related by the simple equation 

Equation: 


Y(z) = H(z)X(z) 


The z-transform decomposes x(7) into its various components along z” which 
passing through the system simply multiplies that value time H(z) and the 
inverse z-transform recombines the components to give the output. This explains 
why the z-transform is such a powerful operation in linear discrete-time system 
theory. Its kernel is the eigenvector of these systems. 


The z-transform of the impulse response of a system is called its transfer function 
(it transfers the input to the output) and multiplying it times the z-transform of the 


input gives the z-transform of the output for any system and signal where there is 
a common region of convergence for the transforms. 


Frequency Response of Discrete-Time Systems 


The frequency response of a Discrete-Time system is something experimentally 
measurable and something that is a complete description of a linear, time- 
invariant system in the same way that the impulse response is. The frequency 
response of a linear, time-invariant system is defined as the magnitude and phase 
of the sinusoidal output of the system with a sinusoidal input. More precisely, if 
Equation: 


x(n) =cos (wn) 


and the output of the system is expressed as 
Equation: 


y(n) = M(w) cos (wn + y(w)) + T(n) 


where T'(7) contains no components at w, then M(w) is called the magnitude 
frequency response and y(w) is called the phase frequency response. If the 
system is causal, linear, time-invariant, and stable, T’(m) will approach zero as 

nm —> oo and the only output will be the pure sinusoid at the same frequency as the 
input. This is because a sinusoid is a special case of z” and, therefore, an 
eigenvector. 


If z is a complex variable of the special form 
Equation: 


z=e” 
then using Euler's relation of e?” =cos (x) + j sin (x), one has 
Equation: 


a(n) = e!”" =cos (wn) + j sin (wn) 


and therefore, the sinusoidal input of (3.22) is simply the real part of z” for a 
particular value of z, and, therefore, the output being sinusoidal is no surprise. 


Fundamental Theorem of Linear, Time-Invariant Systems 


The fundamental theorem of calculus states that an integral defined as an inverse 
derivative and one defined as an area under a curve are the same. The 
fundamental theorem of algebra states that a polynomial given as a sum of 
weighted powers of the independent variable and as a product of first factors of 
the zeros are the same. The fundamental theorem of arithmetic states that an 
integer expressed as a sum of weighted units, tens, hundreds, etc. or as the 
product of its prime factors is the same. 


These fundamental theorems all state equivalences of different ways of 
expressing or calculating something. The fundamental theorem of linear, time- 
invariant systems states calculating the output of a system can be done with the 
impulse response by convolution or with the frequency response (or z-transform) 
with transforms. Stated another way, it says the frequency response can be found 
from directly calculating the output from a sinusoidal input or by evaluating the z- 
transform on the unit circle. 

Equation: 


Zh (n)5| Au) 


z=elw 
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Multidimensional Systems 


Sampling, Up--Sampling, Down--Sampling, and Multi--Rate 


A very important and fundamental operation in discrete-time signal processing is 
that of sampling. Discrete-time signals are often obtained from continuous-time 
signal by simple sampling. This is mathematically modeled as the evaluation of a 
function of a real variable at discrete values of time [link]. Physically, it is a 
more complicated and varied process which might be modeled as convolution of 
the sampled signal by a narrow pulse or an inner product with a basis function or, 
perhaps, by some nonlinear process. 


The sampling of continuous-time signals is reviewed in the recent books by 
Marks [link] which is a bit casual with mathematical details, but gives a good 
overview and list of references. He gives a more advanced treatment in [link]. 
Some of these references are [link], [link], [link], [link], [link], [link], [link]. 
These will discuss the usual sampling theorem but also interpretations and 
extensions such as sampling the value and one derivative at each point, or of non 
uniform sampling. 


Multirate discrete-time systems use sampling and sub sampling for a variety of 
reasons [link], [link]. A very general definition of sampling might be any 
mapping of a signal into a sequence of numbers. It might be the process of 
calculating coefficients of an expansion using inner products. A powerful tool is 
the use of periodically time varying theory, particularly the bifrequency map, 
block formulation, commutators, filter banks, and multidimensional 
formulations. One current interest follows from the study of wavelet basis 
functions. What kind of sampling theory can be developed for signals described 
in terms of wavelets? Some of the literature can be found in [link], [link], [link], 
[link], [link]. 


Another relatively new framework is the idea of tight frames [link], [link], [link]. 
Here signals are expanded in terms of an over determined set of expansion 
functions or vectors. If these expansions are what is called a tight frame, the 
mathematics of calculating the expansion coefficients with inner products works 
just as if the expansion functions were an orthonormal basis set. The redundancy 
of tight frames offers interesting possibilities. One example of a tight frame is an 
over sampled band limited function expansion. 


Fourier Techniques 


We first start with the most basic sampling ideas based on various forms of 
Fourier transforms [link], [link], [link]. 


The Spectrum of a Continuous-Time Signal and the Fourier Transform 


Although in many cases digital signal processing views the signal as simple 
sequence of numbers, here we are going to pose the problem as originating with 
a function of continuous time. The fundamental tool is the classical Fourier 
transform defined by 

Equation: 


Pee i f(t)e Hat 


and its inverse 
Equation: 


oil 
On 


f(t) / F (w) ei dew. 


where 7 = / —1. The Fourier transform of a signal is called its spectrum and it 
is complex valued with a magnitude and phase. 


If the signal is periodic with period f(t) = f(¢ + P), the Fourier transform does 
not exist as a function (it may as a distribution) therefore the spectrum is defined 
as the set of Fourier series coefficients 

Equation: 


P 
ane 5 / f(t) e-2mt/P a 


with the expansion having the form 
Equation: 


f(t) = S_C(k) ePaH/P, 


k 


The functions g;, (t) = eJ2nkt/P form an orthogonal basis for periodic functions 
and [link] is the inner product C (k) = (f (t), gx (t)). 


For the non-periodic case in [link] the spectrum is a function of continuous 
frequency and for the periodic case in [link], the spectrum is a number sequence 
(a function of discrete frequency). 


The Spectrum of a Sampled Signal and the DT FT 


The discrete-time Fourier transform (DTFT) as defined in terms samples of a 
continuous function is 
Equation: 


Fu(w) = S° f (Tn) et 


and its inverse 
Equation: 


can be derived by noting that Fg (w) is periodic with period P = 27/T and, 
therefore, it can be expanded in a Fourier series with [link] resulting from 
calculating the series coefficients using [link]. 


The spectrum of a discrete-time signal is defined as the DTFT of the samples of 
a continuous-time signal given in [link]. Samples of the signal are given by the 
inverse DTFT in [link] but they can also be obtained by directly sampling f(t) in 
[link] giving 

Equation: 


which can be rewritten as an infinite sum of finite integrals in the form 
Equation: 


1 Qn /T 
jiaye - S° / F (w + 2n0/T) efet2r4/T)I duy 
TT ? 0 
Equation: 
1 ae i(w+2ne/T)T 
— wW+2T n 
re Deena) e/ dw 


where F, (w) is a periodic function made up of shifted versions of F'(w) 
(aliased) defined in [link] Because [link] and [link] are equal for all TJ’n and 
because the limits can be shifted by 7/T without changing the equality, the 
integrands are equal and we have 

Equation: 


Pa) = = SF (w+ 2nt/T) = =F, ds), 
L 


where F’, (w) is a periodic function made up of shifted versions of F'(w) as in 
[link]. The spectrum of the samples of f(t) is an aliased version of the spectrum 
of f(t) itself. The closer together the samples are taken, the further apart the 
centers of the aliased spectra are. 


This result is very important in determining the frequency domain effects of 
sampling. It shows what the sampling rate should be and it is the basis for 
deriving the sampling theorem. 


Samples of the Spectrum of a Sampled Signal and the DFT 


Samples of the spectrum can be calculated from a finite number of samples of 
the original continuous-time signal using the DFT. If we let the length of the 
DFT be WN and separation of the samples in the frequency domain be A and 
define the periodic functions 


Equation: 


and 
Equation: 


then from [link] and [link] samples of the DTFT of f(T'n) are 
Equation: 


Equation: 
N-1 . 
=T S° f (In +TNm) eFACn+TNm)k 
m n=0 

Equation: 

N-1 . 

ry i S- f (Tn 4 TNm) ge Nine 

n=0 | m 
therefore, 
Equation: 


F, (Ak) = DFT {fy (Tn)} 


if AT'N = 27. This formula gives a method for approximately calculating 
values of the Fourier transform of a function by taking the DFT (usually with the 
FFT) of samples of the function. This formula can easily be verified by forming 
the Riemann sum to approximate the integrals in [link] or [link]. 


Samples of the DTFT of a Sequence 


If the signal is discrete in origin and is not a sampled function of a continuous 
variable, the DTFT is defined with T’ = 1 as 


Equation: 
Ww) = S° h(n)e 


with an inverse 
Equation: 


1 . 
h(n) = = / H(w) e?*” dw. 


If we want to calculate H(w), we must sample it and that is written as 


Equation: 
Ajay nanje™ 


which after breaking the sum into an infinite sum of length-/V sums as was done 
in [link] becomes 
Equation: 


N- 
H (Ak) =5~ h(n + Nm) e 34% 
m n=0 


M 


if A = 27/N. This allows us to calculate samples of the DTFT by taking the 
DFT of samples of a periodized h(n). 
Equation: 


H (Ak) = DFT {hy (n)}. 


This a combination of the results in [link] and in [link]. 


Fourier Series Coefficients from the DFT 


If the signal to be analyzed is periodic, the Fourier integral in [link] does not 
converge to a function (it may to a distribution). This function is usually 
expanded in a Fourier series to define its spectrum or a frequency description. 
We will sample this function and show how to approximately calculate the 
Fourier series coefficients using the DFT of the samples. 


Consider a periodic signal f (t) = f (t + P) with N samples taken every T 
seconds to give T'n (t) for integer n such that NT’ = P. The Fourier series 
expansion of f (t) is 

Equation: 


ms t) = S- C'(k) e2mkt/P 


k=—0o 


with the coefficients given in [link]. Samples of this are 
Equation: 


_ S- C (k) e2tkTn/P = S- C (k) e2mkn/N 
k=—0o 


k=—0o 


which is broken into a sum of sums as 
Equation: 


co N-1 N-1 
f (Tn) = S- C (k 4 Neg) ermlktNE)n/N _ 2 
£—o0 k=0 k=0 


So C(k+.NQ)| e2lN, 


L-oo 


But the inverse DFT is of the form 
Equation: 


F (Tn) = — SF (k) e®m4IN 


therefore, 
Equation: 


QFI {F(Tn)} =NS~C(k+ Ne) = NO, (h). 
L 


and we have our result of the relation of the Fourier coefficients to the DFT of a 
sampled periodic signal. Once again aliasing is a result of sampling. 


Shannon's Sampling Theorem 


Given a signal modeled as a real (sometimes complex) valued function of a real 
variable (usually time here), we define a bandlimited function as any function 
whose Fourier transform or spectrum is zero outside of some finite domain 
Equation: 


|F(w)| =0 for jw] > W 


for some W < oo. The sampling theorem states that if f(t) is sampled 
Equation: 


fs(n) = f (Tn) 


such that T’ < 22/W, then f(t) can be exactly reconstructed (interpolated) from 
its samples f, (n) using 
Equation: 


{oy Ae 


n=—OCO 


aw 


This is more compactly written by defining the sinc function as 
Equation: 


which gives the sampling formula Equation 53 from Least Squared Error Design 
of FIR Filters the form 
Equation: 


iia S° f; (nm) sinc (xt /T — an). 


The derivation of Equation 53 from Least Squared Error Design of FIR Filters or 
Equation 56 from Least Squared Error Design of FIR Filters can be done a 
number of ways. One of the quickest uses infinite sequences of delta functions 
and will be developed later in these notes. We will use a more direct method now 
to better see the assumptions and restrictions. 


We first note that if f(¢) is bandlimited and if T < 27/W then there is no 
overlap or aliasing in F, (w). In other words, we can write [link] as 
Equation: 


1 oe) 1 n/T 
f (6) / F (w) e3* dw = — F,, (w) e3* dw 
—oo 27 —n/T 


but 
Equation: 


F, (w) = », F (w+ 2rt/T) = ‘>> f (Tn) eT" 


£ 


therefore, 
Equation: 


ro-¢f 


20 —1/T 


T S- f (Tn) asl el" dus 


Equation: 


n —1/T 
Equation: 
sin (2t—an 
=D) 
n T 


which is the sampling theorem. An alternate derivation uses a rectangle function 
and its Fourier transform, the sinc function, together with convolution and 
multiplication. A still shorter derivation uses strings of delta function with 
convolutions and multiplications. This is discussed later in these notes. 


There are several things to notice about this very important result. First, note that 
although f(t) is defined for all t from only its samples, it does require an infinite 
number of them to exactly calculate f(t). Also note that this sum can be thought 
of as an expansion of f(t) in terms of an orthogonal set of basis function which 
are the sinc functions. One can show that the coefficients in this expansion of 

f (€) calculated by an inner product are simply samples of f(t). In other words, 
the sinc functions span the space of bandlimited functions with a very simple 
calculation of the expansion coefficients. One can ask the question of what 
happens if a signal is “under sampled". What happens if the reconstruction 
formula in Equation 12 from Continuous Time Signals is used when there is 
aliasing and Equation 57 from Least Squarred Error Design of FIR Filters is not 
true. We will not pursue that just now. In any case, there are many variations and 
generalizations of this result that are quite interesting and useful. 


Calculation of the Fourier Transform and Fourier Series using the 
FFT 


Most theoretical and mathematical analysis of signals and systems use the 
Fourier series, Fourier transform, Laplace transform, discrete-time Fourier 
transform (DTFT), or the z-transform, however, when we want to actually 
evaluate transforms, we calculate values at sample frequencies. In other words, 
we use the discrete Fourier transform (DFT) and, for efficiency, usually evaluate 
it with the FFT algorithm. An important question is how can we calculate or 
approximately calculate these symbolic formula-based transforms with our 


practical finite numerical tool. It would certainly seem that if we wanted the 
Fourier transform of a signal or function, we could sample the function, take its 
DFT with the FFT, and have some approximation to samples of the desired 
Fourier transform. We saw in the previous section that it is, in fact, possible 
provided some care is taken. 


Summary 


For the signal that is a function of a continuous variable we have 
Equation: 
FT: f(t) + F(w) 
DTFT: f(Tn) - FF, (vw) =D F (wt 20t/T) 
DFT: f,(In) > 7 F,(Ak) for ATN = 2n 


For the signal that is a function of a discrete variable we have 
Equation: 

DTFT: h(n) — Aw) 

DFT: h,(n) — H(Ak)forAN = 27 


For the periodic signal of a continuous variable we have 
Equation: 

C(k) 
NC,(k) for TN = P 


For the sampled bandlimited signal we have 
Equation: 

Since: f(t) — f(Tn) 

f(t) =>¢,, f(Ln) sinc(2at/T — an) 

if F(w) =0 for |w| > 27/T 


These formulas summarize much of the relations of the Fourier transforms of 
sampled signals and how they might be approximately calculate with the FFT. 
We next turn to the use of distributions and strings of delta functions as tool to 
study sampling. 


Sampling Functions — the Shah Function 


Th preceding discussions used traditional Fourier techniques to develop 
sampling tools. If distributions or delta functions are allowed, the Fourier 
transform will exist for a much larger class of signals. One should take care 
when using distributions as if they were functions but it is a very powerful 
extension. 


There are several functions which have equally spaced sequences of impulses 
that can be used as tools in deriving a sampling formula. These are called “pitch 
fork" functions, picket fence functions, comb functions and shah functions. We 
start first with a finite length sequence to be used with the DFT. We define 
Equation: 


where V = LM. 
Equation: 


N=1[ E21 | 
DFT {Uy (n)} = | b(n — um) ra i cd 
Equation: 
Ea Nei | 
= bs d(n — Mine? 


Equation: 


L-1 
e 2aMmk/N -_ 5 eens 
Equation: 


2 L <k>L=0 
~ 10. otherwise 


Equation: 


For the DTFT we have a similar derivation: 


Equation: 
L-1 
DTFT {lly (n 3 soa 5(n — Mm) | aden 
n=—oo |m=0 
Equation: 
—1 foe) 
= S- | s(n Mn) 
m=0 [n=—0o 
Equation: 
L-1 
= e JeMm 
m=0 
Equation: 


_ fL w=k2r/M 
~ 10. otherwise 


Equation: 


M-1 
I=0 


where K is constant. 


An alternate derivation for the DTFT uses the inverse DTFT. 


Equation: 


Equation: 


Equation: 


Equation: 


Equation: 


Therefore, 
Equation: 


1 = 
IDTFT {Uru (w) } = = / Hors (w) e7*”" du 


20 


= | S56 (w — 2nl/M)ei*” dw 
Tv 


= >| 6 (w — 2nl/M) e3" dw 
Tv 


1S * gatn/M ad n=M 
27 otherwise 


I=0 


For regular Fourier transform, we have a string of impulse functions in both the 
time and frequency. This we see from: 


Equation: 
Fr {ur ()}= f oe (t — nT) e dt = > / 3 (t — nT) e ™ dt 
Equation: 
= ai ‘e nae 
Equation: 


2 
= = Dont (~) 


The multiplicative constant is found from knowing the result for a single delta 
function. 


These “shah functions" will be useful in sampling signals in both the continuous 
time and discrete time cases. 


Up-Sampling, Signal Stretching, and Interpolation 


In several situations we would like to increase the data rate of a signal or, to 
increase its length if it has finite length. This may be part of a multi rate system 
or part of an interpolation process. Consider the process of inserting M — 1 
zeros between each sample of a discrete time signal. 

Equation: 


»(n) 7 fone <n>M=0(or n =kM) 


otherwise 


For the finite length sequence case we calculate the DFT of the stretched or up— 
sampled sequence by 
Equation: 


n=0 
Equation: 
MN-1 
Cs(k) = > @(n/M) Un (n) Witw 
n=0 


where the length is now NM andk = 0,1,---, NM — 1. Changing the index 
variable n = Mm gives: 
Equation: 


which says the DFT of the stretched sequence is exactly the same as the DFT of 
the original sequence but over M periods, each of length NV. 


For up—sampling an infinitely long sequence, we calculate the DTFT of the 
modified sequence in Equation 34 from FIR Digital Filters as 
Equation: 


Cs (w) = D7 #(n/M) Un (n) e = fa (m) oem 
Equation: 


= C(Mw) 


where C'(w) is the DTFT of x(n). Here again the transforms of the up—sampled 
signal is the same as the original signal except over M periods. This shows up 
here as C’, (w) being a compressed version of M periods of C(w). 


The z-transform of an up—sampled sequence is simply derived by: 
Equation: 


Equation: 


which is consistent with a complex version of the DTFT in [link]. 


Notice that in all of these cases, there is no loss of information or invertibility. In 
other words, there is no aliasing. 


Down-Sampling, Subsampling, or Decimation 


In this section we consider the sampling problem where, unless there is sufficient 
redundancy, there will be a loss of information caused by removing data in the 
time domain and aliasing in the frequency domain. 


The sampling process or the down sampling process creates a new shorter or 
compressed signal by keeping every M‘" sample of the original sequence. This 
process is best seen as done in two steps. The first is to mask off the terms to be 
removed by setting 1 — 1 terms to zero in each length-M block (multiply x(n) 
by Ij (n)), then that sequence is compressed or shortened by removing the 

M — 1 zeroed terms. 


We will now calculate the length L = N/M DFT of a sequence that was 
obtained by sampling every M terms of an original length-N sequence x(n). We 
will use the orthogonal properties of the basis vectors of the DFT which says: 
Equation: 


Sy ~j2nnl/M M_ if nis an integer multiple of M 
e = 
= 0 otherwise. 


We now calculate the DFT of the down-sampled signal. 
Equation: 


where N = LM andk = 0,1,...,L-1. This is done by masking x(n) . 


Equation: 
N-1 
C4 («| =e (:) tM @ wrk 
n=0 
Equation: 
N-1 , 
= xin ao e J2nnl/M e J2nnk/N 
n=0 M 15 
Equation: 
y MalN=1 
5. x () eJ2m(k+Ll)n/N 
M I=0 n=0 
Equation: 


The compression or removal of the masked terms is achieved in the frequency 
domain by using k = 0,1,...,L-1 This is a length-L DFT of the samples of x(n). 
Unless C'(k) is sufficiently bandlimited, this causes aliasing and x(7) is not 
unrecoverable. 


It is instructive to consider an alternative derivation of the above result. In this 
case we use the IDFT given by 
Equation: 


1 —n 
x(n) = . ce as . 
The sampled signal gives 
Equation: 
, Nal 
y(n) = 2(Mn) = 7 C(k) Ww". 
k=0 


forn = 0,1,---,2— 1. This sum can be broken down by 
Equation: 


, dM a 
C(k+ LI Wyner) 
k=0 1=0 
Equation: 
y tot fMa 
= ye So C(k+L)|Wy™™. 
=0 | 1=0 
From the term in the brackets, we have 
Equation: 
M-1 
Gy(k) = C (k + Ll) 


as was obtained in [link]. 


Now consider still another derivation using shah functions. Let 
Equation: 


xz; (n) = Uy (n) a(n) 


From the convolution property of the DFT we have 


Equation: 
C; (k) = L Uy (k)*C (k) 
therefore 
Equation: 
M-1 
C,(k) = § C(k+ Ll) 


which again is the same as in [link]. 


We now turn to the down sampling of an infinitely long signal which will require 
use of the DTFT of the signals. 


Equation: 
C, (w) = S° a (Mm) em 
Equation: 
— S° a(n) IIy (n)e 7" 
Equation: 
Dew | gd 
= x (n) Se, yan e Jun 
n M 
Equation: 


1 M-1 


= aE S- Soa (n) e Jw 2al/M)n 


l= n 


Equation: 


, Me 
=F C (w — 2nl/M) 


[= 


which shows the aliasing caused by the masking (sampling without 
compression). We now give the effects of compressing x, (n) which is a simple 
scaling of w. This is the inverse of the stretching results in [link]. 

Equation: 


, 
C; (w) = — C (w/M — 2nl/M). 
M > 
In order to see how the various properties of the DFT can be used, consider an 
alternate derivation which uses the IDTFT. 
Equation: 


x(n) = =f CWet de 


TT 


which for the down-sampled signal becomes 
Equation: 


1 7 
xz (Mn) = =| C (w) ee” dus 


TT 


The integral broken into the sum of M sections using a change of variables of 
w = (w; + 2nl)/M giving 
Equation: 


1 Mel ps | 
z (Mn) = = > / C (wy /M + 271 /M) eFi/M+2nl/M)Mn gy), 
=o /—" 


which shows the transform to be the same as given in Equation 9 from 


Still another approach which uses the shah function can be given by 
Equation: 


xz; (n) = Hy (n) a(n) 


which has as a DTFT 


Equation: 
Cx w) = (FF ) Wagar ()*C ) 
Equation: 
sa is C (w + 2nl/M) 
M 


which after compressing becomes 
Equation: 


Approximation Filters. 


Now we consider the effects of down—sampling on the z-transform of a signal. 
Equation: 


Applying this to the sampled signal gives 
Equation: 


A;(Z)= Soa (Mn) 2!" = SE: (n) Hy (n)z ” 


Equation: 


M-1 — ) 
i x (n) elemnl M gt 
Equation: 
M-1 a 
= Vem ferr 2} 
I=0 7 
Equation: 
M-1 
= x (ee z) 


which becomes after compressing 
Equation: 


= 
= x (eon zt), 
= 


This concludes our investigations of the effects of down—sampling a discrete— 
time signal and we discover much the same aliasing properties as in sampling a 
continuous-time signal. We also saw some of the mathematical steps used in the 
development. 


More Later 


We will later develop relations of sampling to multirate systems, periodically 
time varying systems, and block processing. This should be a very effective 
formulation for teaching as well as research on these topics. 


